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Abstract 

The plant hormone auxin is fundamental for plant growth, and its spatial distribution in plant 
tissues is critical for plant morphogenesis. We consider a leading model of the polar auxin flux, and 
study in full detail the stability of the possible equilibrium configurations. We show that the critical 
states of the auxin transport process are composed of basic building blocks, which are isolated in a 
background of auxin depleted cells, and are not geometrically regular in general. The same model 
was considered recently through a continuous limit and a coupling to the von Karman equations, to 
model the interplay of biochemistry and mechanics during plant growth. Our conclusions might be 
of interest in this setting, since, for example, we establish the existence of Lyapunov functions for the 
auxin flux, proving in this way the convergence of pure transport processes toward the set of critical 
configurations. 



1 Introduction 



The plant hormone auxin plays a fundamental role in plant development ( Reinhardt et al.^|2000 20031, 
and its spatial distribution in plants tissues is critical for plant morphogenesis. Auxin accumulation is 



spatially localized in specific set of cells, where it induces the emergence of new primordia (Reinhardt 
et al.,||2000||. A fundamental problem consists in understanding how such auxin maxima appear, and 

pOOQl). On the 



Hamant and Traas. 



how they induce the regular pattern observed in plants (see e.g. 
other hand, experiments show that phyllotaxis strongly depends on the plant physical properties, more 
precisely on elasticity ( Green, [1980 Dumais and Steele, |2000[ Dumais,|2007 1, and physical forces provide 
information for plant patterning (Hamant and Traas, 20091. Basically, turgor pressure induces stress, 



which is related to the associated deformation or strain through Young constants: see e.g. |Boudaoud,| 
(20101 where these notions are explained in the context of plant growth. Experiments have shown that 



lowering the stiffness of cell walls in the meristem leads to the emergence of new primordia (Hamant 
et al.,|p008 l. However, the interactions between physics-based and biochemical control of phyllotaxis is 
still poorly understood. 



Recently, new biologically plausible mathematical models of auxin transport have been proposed (Bar 



bier de Reuille et al.,|2006 


Heisler,|2006 


Jonsson et al.,|2006 


Smith et al.,|2006 



to reproduce some aspects of phyllotaxis in simulations. New mathematical models were also proposed for 



2008) 



the plant mechanics (Mjolsness, 20061, and for the interaction between mechanics and biochemistry (Ship- 
man and Newell,]|2005 Newell et al.. 



In the latter, the authors use the model for the polar auxin 
flux proposed in I Jonsson et al., ( 2006 1 for modelling the stress field in their mechanical model. It should 
be stressed that all these models are based on hypotheses that have not been verified experimentally, 
however they provide new scenari for understanding plant growth that can be tested experimentally. 
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Figure 1: Inflorescence shoot apical meristem of Arabidopsis thaliana. Zones with high auxin concen- 
tration are highlighted by the fluorescent yellow signal auxin reporter from DR5::YFP. The red signal is 
highlighting cell walls stain, using propidium iodide. 



Auxin occurs in various plant tissues, where it is transported by polar cellular transport in various 



et al., 


(2006 


1 or 


Bayer et al., 


(2009 



In the following, we consider the models in Jonsson et al., (2006) and Smith et al., (20061, based on 
polar auxin flux. Polar auxin flux results from uneven accumulation of the auxin transport regulator 
PIN in cell membranes. An essential component is a positive feedback between auxin flux and PIN 
localization, resulting in the reinforcement of polar auxin transport to dedicated routes which develop 
into vascular tissues. We will not enter here into these considerations, but focus on simple models of 



transport processes (see e.g. the discussion in Jonsson et al., (2006) and Shipman and Newell, ( 2005[ )), 
where a quasi-equilibrium is assumed for PIN proteins. The molecules present in some cell i may be 
transported to any neighbouring cell j, but they are preferentially transported to the neighbours with 
the highest auxin concentrations. 

Traditionally, models of patterning and morphogenesis have used reaction-diffusion theory. Turing demon- 



strated how, under some hypotheses, the regular patterns observed in phyllotaxis can be predicted (Tur- 
g^|1952 1 . He showed that a combination of diffusion and a chemical reaction could give rise to regular 



patterns. Interesting models are described in Meinhardt, (19821; Thornley, (1975) which can, under some 
hypotheses, predict phyllotactic patterns. As stated previously, the auxin flux is strongly polarized, a 
phenomenon that cannot be described with reaction-diffusion models. The recent mathematical models 
given in 



Barbier de Reuille et al., (20061; Jonsson et al., (20061; Smith et al., (2006) are based on trans 



port processes. Mathematically, mass transport processes are not well understood, and their study is a 
challenging problem. We propose here a mathematical study of related dynamical systems. We focus on 
their critical points and analyse their geometrical structure and stability. 

Besides stable auxin peaks, the model generates intervening areas of auxin depletion, as it is observed 
experimentally. These auxin depleted sites reflect an indirect repulsion mechanism since auxin molecules 
diffusing through the tissue will be attracted to the peaks, and diverted from the depleted areas. This 
idea of repulsion or spacing mechanism was already considered a long time ago ( Hofmeister,||1868 ) . 

The auxin flux is present everywhere in the plant, so that we choose to describe the various plant cells 
as a connected graph (A, The node set A represents the cells and E the set of edges. Any edge 
e = (i — > j), i, j € A indicates that some auxin molecule can move from cell i to cell j. This graph is 
undirected, and we write i ^ j to denote that cells i and j are nearest neighbours, so that auxin can move 
from cell i to cell j , at some rate qij . These transition rates are not well understood at present time and 
one must rely on simple models. They should capture the fact that an auxin molecule present in some 
cell i has the tendency to move to a cell j ^ i when the concentration aj of auxin molecules present in 
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cell j is high. The simplest model accounting for this idea is given by 



Jonsson et al., 



(20061 



ak 



(1) 



for some positive constant k, which is of Michaelis-Menten or Monod type. Let L = |A| be the number 
of cells. In the model given in Jonsson et al., ( 2006[ ) (see also Smith et al., (20061; Sahlin et al., (2009l), 
ai{t), for i — 1, • • • ,L, denotes the concentration or the number of auxin molecules in cell i at time t, 
and is assumed to evolve according to the differential equations 



^Ma) = Dj2iak-ai)+Tj2 



ak 



ak 



^qki{a) 



=<?ifc(a) 



(2) 



for i = 1, . . . , L. The term aiQik gives the mean number of auxin molecules moving from cell i to cell k 
per unit time, and D Efe~i('^fc ~ ^ diffusive part, usually assumed to be weak with a small diffusion 
coefficient D. The second term corresponds to the mass transport process, which is known to be the 
main actor of the patterning process in plants. One can add auxin production and degradation terms, 
but, there is no clear biological evidence about where auxin is produced, and experiments show that it is 
not produced in the meristem, but imported from the leaves (Reinhardt et al.^|2000 20051. 



1.1 Results 

Direct quantitative measurements of auxin distribution in plant tissues are very difficult due to the small 
size of the meristematic tissues at the time of patterning. Therefore, biologists rely on indirect markers 
based on auxin-regulated genes that encode fluorescent proteins. Figure [T] shows a typical output, where 
domains rich in auxin appear as regions of strong green fluorescence. The pattern is quite noisy; this 
might be due either to the indirect experiments, or to the fact that the number of auxin molecules is 
not too high. ([2| might model the limiting behavior of this random particle system when the number of 
molecules tends to infinity. We introduce such a particle system in Section |2] and justify equations like 
(|2| using law of large numbers. 



We then focus on the properties of ([2|, like the non-negativity of the solutions (see Proposition 3.1 1. This 
dynamical system can be written in the compact form 

da 

dt = 

where a{t) — (ai(i))i<i<L is the vector of auxin concentrations. The related critical points are the vectors 
a* satisfying f{a*) = 0. They are the candidates for describing the equilibrium auxin concentrations. 
For example, ai — Q means that there is (almost) no auxin molecules in cell i, while a subset of cells J 
such that aj > for j £ J indicates a hot spot which might correspond to an auxin peak. 

The critical points play a fundamental role in the dynamic, and one can suspect that any solution a{t) of 
([2]) will approach such critical points as t is large. Of course, this is wrong for general dynamical systems, 
but here, the model is supposed to catch pieces of biological reality, and the robustness of the regular 
geometries observed in plants suggests that this might well be the case. Some of these critical points are 
repulsive or unstable, that is, the orbits or the solutions of (|2| will avoid them. In the contrary, some of 
them will be attractive. Given a critical point a*, a mathematical way of checking the stability or the 
unstability of a* is to compute the Jacobian d/(a*), by retaining only its spectrum, that is the set of 
all eigenvalues of d/(a*). For example, a* is unstable when there is an eigenvalue having a positive real 
part. 

Definition 1.1 We say that a critical point a is stable when all the eigenvalues of the Jacobian evaluated 
at a have non-positive real parts. 
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Section [5] is concerned with the characterization of the set of critical points, mainly focusing on pure 
transport processes. 



For D = 0, we first consider critical points a > 0, meaning that > for all i. Corollary |5.3| shows that 
such elements are precisely the positive solutions of the linear equation 



Ta = c 1, c constant, 



(3) 



where F is the adjacency matrix of the graph G, with entries F^j e {0, 1} such that Fy = 1 if and only 
if cells i and j are nearest neighbours, and 1 is the vector having all components equal to 1. 

Next, we focus on critical points such that Ui — for i belonging to some subset / C A = {I,-- - ,L}. 
They correspond to auxin depleted cells. The graph decomposes into a product of sub-graphs 7, which 
are the connected components of the sub-graph of G induced by the node set J = A\I. We thus look for 
a having positive components aj > for j S J, which should correspond in some sense to auxin peaks. 
We obtain the distribution of auxin in such components, denoted by a|^, by solving the linear systems 
F^a|^ = c^l|^. A typical example of such configurations is given in Figure [s] where the elements of I are 
black and the various components 7 red. 

We then turn to the asymptotic behavior of the solutions of system and establish in Proposition 
|6.3| that every solution converges toward the set of critical points. Our technique is based on Lyapunov 
functions, that is, we look for a function H{a) which should be decreasing along the orbits of Q, like 
energy in physics. We proved that, for pure transport processes with D = 0, the function 

H{a) = -K{l,a) - i(a,Fa), 



where (•, •) denotes the scalar product, satisfies 



dH{a{t)) 
dt 



< 



for any solution of ([2| . [Newell et al., ( 2008 1 also considered the differential system Q by taking a spatial 
continuous limit, and showed that the limiting equation is a p.d.e. similar to the von Karman equations 
from nonlinear elasticity theory: 



dw 

Ik 



PAw + const 



nonlinear terms. 



The von Karman equations are of gradient type (see e.g. Shipman and Newell 



potential is given by the elastic energy. These energy functionals were then used in Newell et al. 



(|2005[)), where the 
p08l 



and Newell and Shipman, ( 2005 ) to provide a very interesting mechanical explanation of the appearance 

However, the limiting equations associated 

For the basic 



of Fibonacci numbers in plant patterns based on buckling. 

with the auxin flux are not of gradient type, see the discussion in 'Newell et al.,| (12008 1 
dynamical system ([2|, our result shows that the system is minimizing the energy H, without being of 
gradient type. 

Section [7] considers stability, and Proposition 
at a is given by 



7.1 



shows that the Jacobian d/(a) = (dfi/daj) evaluated 



d/(a) = 



1 

7V2 



2d(a)F(c id-d(a)F 



where d{a) is the diagonal matrix of diagonal given by a. This permits to check the stability of the 
critical points for various graphs. We present various results on graphs of interest for plant patterning 
questions, like the circle or the two-dimensional grid. As stated previously, the positive solutions a\^ to 
the linear system F^a|^ = c-yl|^ provide restrictions of the critical points to the connected components 
7. We give a particularly simple condition on the sub-graph 7 of G induced by the set J = A\I ensuring 
the non-stability of a|-y. Let Mi, i G A be the neighbourhood of i, that is the set of nodes j such that 
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j i and j ^ i. The configuration a\.y is unstable when the sub-graph 7 contains a path of length 4, of 
the form 

io ii 12 ^ h, 



such that 



ii G 12 ^ and ig ^ TV,, 



4" n 



Figure 2: Example of components 7 of the two-dimensional grid that can potentially yield stable config- 
urations, see Corollary |7.3| 



For example, if G is a two-dimensional grid, any stable configuration is composed of patches of the basic 
building blocks given in Figure [2j These patterns are however not geometrically regular in general, see 
f igure [3| The more involved mode l of [Smith et al.,| (|2006[ ) , which uses PIN proteins in a direct way (here 
we assume a quasi-equilibrium, see|Jonsson et al7^ (2006 1), produces more regular patterns in simulations. 



In this setting, the transition rates are forced to follow exponential distributions. Hence, a strong selection 
based on rates of the form exp(&ai), 5 > instead of the linear function ai seems to regularize the critical 
points. Of course, it might be interesting to justify such a choice biologically. We also argue in what 
follows that the critical configurations produced by the auxin fiux might be more regular when coupled 
to periodic potentials. 



Figure 3: A potentially stable configuration when the graph G is a rectangular grid, for the pure transport 
process. The black circles correspond to the values a.^ = 0, i G / (the set auxin depleted cells), while the 
red circles are such that > 0, corresponding to auxin peaks. One can construct the set of all stable 
configurations by playing with the building block given by the square, the star, and the various parts of 
the star. This shows that dynamical system ^ does not necessarily produce regular patterns. We can 
however give examples where such configurations are unstable, see Section [73] 



It might well be that the auxin flux self-organize in regular patterns when coupled to mechanical forces, 
for example, as already stated in the Introduction, see Newell et al., (20081. In the same spirit, we 



introduce a simple model coupling the auxin flux to a potential </>, which might model deformations. 
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(a) a{t) for f « (b) a{t) for large t. 

Figure 4: Simulation of the orbits of the differential equation Q with T = 1, D = and a potential 
(j){x, y) = sin(47ra;/j4) sin(47ry/_B) on a torus, where x = 1, • • • , A and y = 1, • • • , -B. The initial state is 
flat, (b) shows the state a{t) for large t: one sees regularly spaced auxin peaks, which are isolated in 
a background of auxin depleted cells. The potential and transport terms drift thus the process toward 
more regular patterns, while the transport process creates domains of auxin depletion. 




curvature or effects related to the meristem elasticity. We provide an example of the form 



dt 



(4) 



1 = 1,. 



at al. 



L. If the potential itself has some regularities, as it is the case in specific model given in |NewelT| 
(20081, the auxin flux will exhibit much more regular patterns, see e.g. Figure [4] 



Finally, the model provides an interesting conclusion: for most graphs, stable configuration are composed 
of building blocks isolated in a sea of auxin depleted cells. This might be the basis for repulsion between 
primordial auxin molecules will not have the tendency to move toward them, leading to indirect repulsion. 



The idea of such repulsive force appeared a long time ago in the work of Hofmeister, ( 1868 1. Many authors 
have used this hypothesis to develop very interesting mathematical models, all leading to phyllotactic 
patterns observed in nature, like Fibonacci numbers, the Golden Angle or helical lattices, see |Adler,| 



(1974 


); 


Atela et al.. 


(2002 


I; 


Douady and Couder, 


(1996 


); Kunz, 


(1995 


); 


Levitov, 


(1991 



2 A stochastic model of auxin transport 

We consider a stochastic process related to differential equation ([2]), describing the random numbers of 
auxin molecules rjt{i) G N present in cell i at time t, i = 1, • • • , The state space of this stochastic 
process is denoted by £7^ = N'^, where A is the set of L cells (the nodes of the graph). Looking at 
equation (|2|, we define transitions by supposing that any auxin molecule present in cell i at time t can 
be transported to a neighboring cell j at rate qij {rf) of the form 

when ri{i) > 1. This defines a Markov process with state space $7^, describing the stochastic moves of 
the various auxin molecules. Let M denote the total number of molecules. It turns out that the ordinary 
differential equation (|2| describes the large M limit of the stochastic process (weak noise limit). This 
random particle system is then described as a gaussian process XM{t) ~ rjt/M in drifted by the 
solution a{t) of (|2| for some covariance function. This approximation will be mathematically rigourous 
if the constants k and R are related in such a way that R = Mk, and the limiting behavior of the rescaled 
number of auxin molecules is such that r]t{i)/M w ai(t), where a(t) solves ([2]), with X]ieA'^»(^) = 
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Such stochastic particle systems are known as density dependent population processes, and the above 
limit has been treated in detail in|Ethier and Kurtz, (1986), and corresponds to a law of large numbers. 
Notice that different kinds of limits can also be considered. Stochastic mass transport processes of this 
type have also appeared in physics, and are known as generalized zero range processes, see e.g. [Evans | 



et al., (20061; Luck and Godreche, (20071; Grosskinsky et al., (20111; Kipnis and Landim, (19991. In this 



setting, hydrodynamical limits are considered, when both M and L tend simultaneously to oo in such a 
way that M — pL, for a fixed density. Simulations show the appearance of condensates when p is larger 
than a critical threshold pc, which might represent auxin peaks in some way. Mathematically, the theory 
of condensation is not developed at present time for these general processes, so that we here focus on the 
weak noise limit. 



The gaussian approximation of rjt/M is defined as follows: for i — \, ■ ■ ■ , L, consider the unit vectors 
with ei{j) — when j ^ i and ei{i) — 1. Let d/ be the Jacobian d/ = {dfi/daj))ij=i^... ^l- 

For simplicity, we illustrate the transition rates for cells arranged along a circle: the rate functions are 
given by functions /3i{a), I ^1} , satisfying 

/3e.^,_e,(a)^(Ja,+r , "'+^"' ) ; for z = 1, . . . , L, 

/3e._i-e.(a) = (Z?a,+r^ ) ; for z = 1, . . . , i, 

/3;(a)=G for ; ; ^ ei_i - Ci, Ci+i - Ci. 

For example e^+i — means that an auxin molecule of cell i has been transported in cell i + 1. For 
arbitrary graphs, the definitions of the rates Pi are similar. 

With these notations, we can define the matrix G 

G{a)= Y,Pi{a)ll*, 

which will be an essential element of the covariance matrix associated with the gaussian approximation. 
Consider the following matrix valued differential equation 

^^M^d/(a(t))0(t,s), 0(s,s) = id. 
Then, as M is large, one gets that (see e.g. ( Ethier and Kurtz"^|1986 1) 



where Vi is a gaussian process of mean O)y(O) and of covariance function 

/•minli.r} 

Cov(t/(t), V{r)) = / <i>(t, 5)G(a(s))0(r, s)*ds. 

Jo 

3 Basic properties of the auxin flux 

Proposition 3.1 Every solution a of ^ starting in Rk.^ remains non-negative, and is conservative, 
that is, 

L L 

\ft e M>o, J2 = = P^- 

i i 

Moreover, the system admits a unique solution defined over [0, +oo). When 0^(0) > 0, then ai{t) > 0, 
Vt > 0. For pure transport processes with D — 0, ai{0) ~ =^ ai{t) = 0, V< > 0. 
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The proof of proposition |3.1| is given in Section |8] 
Let us rewrite the system Q, for I < i < L 

a. = Dj2ak + TY,{ ^ y - °" a " (6) 

with the initial condition a(0) G M:^. 

Proposition 3.2 // the graph is connected and D > 0, the only critical point of ^ in admitting 
zero components is the origin. 



Proof : 

Let a; — where is the i-th component of a critical point a G M;^ of ([6]). Clearly ([6| entails a/c = 

and the non-negativity of each term. Ofe = for all /c i. Since the graph is connected we deduce that 
flfc = for all 1 < fc < L. □ 

Proposition 3.3 Let us assume that the graph is connected and D > 0. If^^^i afc(O) > 0, then for all 
j G {1, L}, we have lim ^^^^a,;(^:) > 0. 



To prove the previous proposition, we will use the following Proposition, see Gabriel et al., (19891. 

Proposition 3.4 Let f : M_|_ — > M 6e twice differentiahle and bounded together with f . If, as n ^ +oo, 
tn t +00 and f{t,,) \jmt^+ocf(t) (or f{t„) ^ \hRt^+oof{t)), then f{t„) ^ 0. 

Remark 3.5 (1) The boundedness of f and f implies the one of f. 

(2) The assumptions in the preceding proposition can be weakened without changing essentially the proof: 
"/ : M-i- — 5- M 6e twice differentiable and bounded together with f " can be replaced by "f : M_|_ — > M 
is bounded and differentiable and f is uniformly continuous ". 



Proof : 

If afc(O) = for all 1 < fc < L, then the unique solution is identically zero. Otherwise J2k=i '^k{0) > 0. 
Let us suppose that for some i G {1, L}, 

Let us introduce the notation — lim ^_,_^^^a, (Q. Since ai{t) is bounded together with its second 
derivative, the preceding proposition applies and for any sequence t„ f +oo such that ai(t„) Qn, we 
have ciiifn) — >■ as n — >■ -|-c». Every ak{tn) being bounded in the right-hand member of the equation 
for ai{tn), we conclude that lim„_i.+oo ^ X]fc~i ~ 0. The non-negativity of each ak{tn) entails 
lim„_>.+oo o,k(tn) — — (kk for every k ^ i. According to the above proposition, lim„_5.+oo OLk{tn) ~ for 
every A; ~ i and since the graph is connected, repeating the same argument provides lim„_i.-j,oo o,j{tn) — 
for every j G {1, L}. Thus = lim„_>._|_oo X]i<j<L '^ji^n) =^ ctfc(O) > 0, a contradiction. 

As a consequence, for D > 0, it is impossible to have limt^+co ai(t) — 0, and thus none of the com- 
partments can become empty asymptotically. □ 
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4 Tools from Markov Chain theory 

We will use notions from Markov chain theory, and hence consider generators Q : A x A — > M, 
Q = {qij, i, j G A}, such that 

Qij > 0, for i ^ j and = - ^ Qt]- 

For example, the auxin flux described by ([2| contains implicitly a generator Q{D, T, a) given by 

qij{D,T,a) ^ D + Tqij{a), i ^ j, 

< q,j{D,T,a)=0, iooj,i^j. (7) 

q,i{D,T,a) = - Ej^j 9u (-D, T, a), 

where we set 

/ ^ 

Q is irreducible when for any pair of nodes there is a path io = i ^ ii ^ 12 ^ ■ ■ ■ ^ ik = j such 

that (?i„i„^i > 0, n = 0, • • • , fc — 1. When Q is irreducible, one can prove that there is a unique invariant 
probability measure tt satisfying 'n*Q — Q. 

An irreducible transition kernel Q of invariant probability measure tt is said to be reversible when 

TTj^y = -Kjqjt, yi ^ j. 



5 Characterization of the critical points 

We can write ^ in the more compact form 

^ =/,(a) = ^(ajq,,(-D,r,a)-a,(j,,(-D,r,a)), ^ = /(a) = a*Q(A T, a). 

Our first aim is to look for the critical points of the above dynamical system, that is, to find the element 
a G solving the equations /(a) = 0, which can be rewritten as a*Q{D, T, a) = 0. Hence, any solution 
to /(a) = is an invariant measure associated with the transition function Q{D,T, a). We will use the 
following facts: 

• When D > 0, the generator Q{D,T,a) is irreducible. 

• For pure transport processes where D = and T > 0, Q(0, T, a) is irreducible if and only if > 
Vi. 

In the irreducible case, let 7r(a) denote the associated positive invariant probability measure. We thus 
look for a > such that 

^—^nia). (8) 

5.1 The irreducible case 
Pure transport processes 

If (5(0, r, a) is reversible, the equation /(a) = is equivalent to the set of equations 

aiqij{0,T,a) = ajqji{0,T,a), i ^ j. (9) 
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In what follows, we will use the functions 



Nk = Nk{a) = K + ^ aj. 



(10) 



Lemma 5.1 Let G he a connected graph. Assume that D ~ and T > 0. Then Q{0,T,a) is reversible 
Va > 0, of invariant probability measure given by 



VZ a / jsA 



where 

Z[a) = QiNi = CLi + aiQj. 
In this case, a > is a critical point with f{a) ^ if and only if Ni{a) does not depend on i, with 



(11) 



E 



(12) 



ieA 



Remark 5.2 The transition rates qij{a) are similar to the rates associated with a family of Markov chains 
used in the study of vertex-reinforced random walks, see Benaim, (1991); Benaim and Tarres, (2008) and 
Pemantle, (1992), and Lemma \5.1\ is an adaptation of these results. Interestingly, such vertex-reinforced 



random walks are approximated by deterministic dynamical systems called replicator dynamics, of the 
form 

^=a,(iV;(a)-i/'(a)), 

where N[{a) = Ni{a) — k and H' {a) — EieA^i-^i- ^^^'^ setting, the function H' plays the role of i 
Lyapunov function. We will also find a similar Lyapunov function, see Section^ 

Proof : 

Assume, without loss of generality, that T = 1. First notice that 

i:.(a),„.(o,r,a) . 1:1^9- 



z{ 



1 \ ^ a, i \ ^ 



ai{Ni - k) 



The identity 



shows that 



7r(a)jfe(0,T, a) 



(«) ' Z{a) 



ai{Ni - k) 



Z(a) ' 

^7r(a)jgj-,(0,T,a) +7r(a),qi,(0,T,a) ^ 0, 



so that 7r(a) is a invariant probability measure for (5(0, T, a), a > is a critical point with f(a) = if 
and only if ^ ° ^ is an invariant measure for Q(0, T, a). Because of the unicity of the invariant measure, 
we obtain 



E 



ieA ' 



□ 



Let r be the adjacency matrix of the graph G = (A, E), that is, the matrix with entries given by V^j = 1, 
when i ^ j and i ^ j, and — otherwise. We summarize the above results in the following 
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Corollary 5.3 (Pure Transport Processes) Assume that D = and T > (no diffusion), and con- 
sider only positive a > 0. Then, 



/(«) 



where 



if and only ifVa — c(a)l, 1 = (1, ■ 
'E^eA Ej~» (a, Ta) 



c(a) 



(a,l> 



(13) 



(14) 



Remark 5.4 Let c be a constant, and let a (if it exists) be such that Ta = cl and a > 0. Then a. is a 
critical point and c is given by {H). 

Example 5.5 (The one-dimensional cycle) Assume that the L cells are arranged on a cycle. The 
pure transport process (D — 0) is reversible, so that the critical points a > of dynamical system ^ 
are solutions of linear system (13). We illustrate some results given in Section 7.4 When L > A is a 
multiple of 4, the set of critical points a G M.^ forms a two dimensional sub-manifold Mc o/M^ given by, 
when p = 1/L, 

Mc = {(ai, 02, -ai + 2/9, -02 + 2/?, ai, 02, -ai + 2p, -02 + 2p, • • • ); at e (0, 2p), k = I, 2}. 

When L > A is not a multiple of 4, Mc is reduced to the uniform configuration Mc = {(p, p, • • • , p)}. We 
will see that the uniform configuration is always unstable, and that the other critical points are unstable 
when a > 0. However, the boundary points are all stable. 



General transport processes 

Lemma 5.6 Assume that G is connected, and that both D and T are positive. For a > 0, /(a) — if 
and only if there exists a constant c such that a solves the following system of quadratic equations: 

{a, - ^)N,{a) + a, = c a,N,{a), i = 1, • • • , A. (15) 

Proof : 

Let fii — {ai — D/T)Ni, i ~ 1, - ■ ■ , A. Then fi = (pi)i<i<A behave 



{lJ.Q{0,T,a))i = ^p,jqji{a) + p,qi,{a) 



which gives the diffusion term contained in /. Hence, one can rewrite the equation /(a) = as 

{fi + a)Q{0,T, a) = 0. 

By assumption, a > so that Q{0,T,a) is irreducible as a Marfcov generator, and hence has only one 
invariant probability measure. The linear space composed of invariant measures is one-dimensional, so 
that the measure /x + a is proportional to 7r(a). The result is a consequence of expression for n{a) given 
in □ 

The next paragraph generalizes the diffusive part to model the effect of potentials on the auxin flux. 
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Inclusion of potentials 



As stated in the Introduction, 
play a role in plant patterning. 



experiments have shown that both mechanical and biochemical processes 



We here adapt some ideas of Newell et al., (^2008 1 and Newell and Shipman, 



(2005) to our discrete setting. The former considered the discrete model (|2| by taking a continuous limit, 



|2008 
jllby 



resulting in a p.d.e. describing the time evolution of auxin concentrations, which is coupled to the von 
Karman equations from elasticity theory. These equations describe the deformations of an elastic shell or 
plate subject to various loading conditions. Usually, the in-plane stress is described using Airy functions 
which are potential for the stress field. Here, we will simply suppose that this potential is given by some 
function {4>i)i<i<L- We also suppose that the auxin flux is directed in part by these potentials and assume 
a model of the form 

dai{t) . . ■r-^ , 



(16) 



i = 1, . . . ,L. We will see in the sequel that the critical points associated to (|2| exhibit regular geometrical 
patterns locally, but not necessarily globally. The potential might be defined in such a way to reproduce 



the patterns obtained when considering mechanical buckling, and the model defined by ( 16 1 might then 
lead to more regularly spaced auxin peaks, see Figure |4] 



Lemma 5.7 Assume 



model of the form (16), with D > and T > 0. Let a > 0. Then fi{a) + 

such that 



'^j^i{cLj4'i — ciiff'j) — ^ <^^d only if there exists a constant c G 

n 1 

i = c aiNi{a), i = 1, 



The proof of Lemma |5.7| is identical to the proof of Lemma |5.6 



5.2 The reducible case 

We can adapt the previous notions to the case D = and reducible transition kernel Q(0,T, a), that is 
when some vanish. In this case, there is a pair of nodes i and j such that 

m 

fe=i 

for all paths 7 : iq = i — )• ii — )• • • • — > im = j taking i to j in the graph G — (A, E). 

Example 5.5 shows that the critical points associated with ^ on a circle form a manifold when L is a 



multiple of 4. We also assert that the boundary points obtained from Mc by setting ai = are stable 
We will thus consider subsets I C {1, ■ ■ ■ , L} corresponding to the sites i where = 0. We will denote 
by a\i the restriction of any a to I. The same notations apply for generators and adjacency matrices 
where one conserves only the transitions rates qij{a) such that i, j G A \ I. According to Lemma 5.1 
these sub-transition kernels are reversible for a such that a\\\j > 0. If one removes the nodes i G I, the 
graphs decomposes as a product of connected components 7, which form the sub-graph of G induced by 
the nodes of J = A \ /. The special form of the vector field associated with ^ ensures however that 
the set of critical values such that a\j = 0, / C {1, • • • , L}, can be obtained by considering a family of 
transitions functions (5^(0, T, a|^). For each component 7, Corollary 5.3 shows that the related critical 
points are obtained by solving linear systems of the form 

r^a|-y = c^l|7, (17) 

where is the adjacency matrix of the sub-graph 7, and the are normalization constants chosen in 
such a way that a.^ = pL. The set of critical points is then obtained by taking the direct product of 
the sets of critical values associated with the sub-graphs 7. 
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6 Asymptotic properties of the auxin flux for pure transport pro- 
cesses 



We consider the convergence of the dynamical system Q when D — using the method of Lyapunov 
functions. Suppose without loss of generality that T = 1. We look for a function H{a) such that 

'i^('^(^».(Vi7(a(t)),^)<0, Vt>0. 



dt 



dt 



If furthermore this function is bounded, then H{a{t)) converges, and we can in this way get useful 
information concerning the convergence (e.g. toward the set of critical points) of a{t) solution of Q. 



Lemma 6.1 Assume that D = and set T ~ 1. Let 

H{a) = ^afc(7Vfc(a) + k) = -K^a^ - \ ^^a]ak, 



(18) 



fceA 



fceA 



fc j~fc 



where the functions Nk{a) have been defined in (10). Let a{t) he a solution of the o.d.e. ^ such that 
ai(0) > 0. Then 

"^"^ = - J E E m^^^^ - N,f < 0, vt > 0. 



(19) 



Notice that 



keAj~k 



dH 



since the function Nk — Nk{a) — n + X)j~fc does not depend on the variable a^. 

Proof : 

One can write 



dH{a{t)) 
dt 



fceA j~fe 



fceA j~fc 



aj ak 



(^fc - N,) 



keA j^k 
1 

2 



fceA:/~fe 



keKj^k 



By Proposition 3.1 0^(0) > 0, Vi, implies that ai[t) > 0, Vi, Vt > 0, so that qkj > and qjk > 0, Vfc ^ j, 
and yt > 0, proving the assertion. □ 



To prove the convergence of the auxin flux, we use a Theorem of Lyapunov- LaSalle (see LaSalle, (1976l). 
Introduce the notation 



BH 1 

^(^) = E = - 9 E E MkiNk - 



fceA j~fc 



Consider the sets 



n = {x€ [0, 2p]^ I ^ = pL} and = {cc G | ij(a;) = 0}. 
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Lemma 6.2 The set is the set of critical points. 



Proof : 

Let X E fl. Then H{x) = if and only if for all pairs j ^ k, either Xj = 0, Xk = or Nj = Nk- Let 
Ix ■— {« G A; Xi — 0}. Then H{x) = if and only if, for all pairs of neighbours j ^ k such that j £ A\Ix 
and k G A\ I^, one has that Nj = Nk- Let 7 be the connected component of the graph containing 



this pair (see Section 5.2 1, with Nj ~ Nk — c^, for some positive constant c^. Then, Ni = c^, Vi G 7. 
One then gets that H{x) = if and only if the function N is constant on the connected components 
7 associated with I^- Hence, for each such component, one has that F The results is a 

consequence of Corollary |5.3| and of the results of Section |5.2[ □ 

Let be the largest invariant subset of E^. As E^ contains only the critical points of /, Eq is invariant. 
Hence, Mn = En- 

Proposition 6.3 Let a{t) be the unique solution of the o.d.e. ^ with a{0) G 51. Then a{t) G 51, Vi > 
and a{t) converges to Mq as t 00. 

Proof : 

Proposition |3 . 1 1 shows that the compact set 51 is invariant. The contin uously di f Ferent iable function H is 



such that H{x) < 0, Va; G 51. The results then follows from a result of LaSalle, (1976). □ 



Corollary 6.4 Every limit point of a trajectory a(t) is a critical point i.e. if for <„ 00, a(t„) Oo 
then Ooo G Mo. 

Proof : 

If Qoo ^ 51, as Eq = Mil is a closed set then d() > 0. It's a contradiction with the proposition 6.3 □ 



Remark 6.5 (Global minimizers of H) The literature contains results on the set fJ,{G) of minimizers 
of H when X^igA '^i — ^- authors of (Motzkin and Straus, ^^1965 ) proved that maxa(a. Fa) = {uj{G) - 



l)/ijj{G), where lo{G) is the clique number of G, that is the order of the largest complete sub-graph of 
G. Moreover, they obtained that the absolute minimum of H is achieved at an interior point of the unit 



simplex if and only if G is a complete multipartite graph. Various results were then obtained in (Waller, 



1971). where for example it is proved that I-l{G) is a simplicial complex, having an automorphism group 



similar to that ofG. In some sense, /i(G) mirrors some of the geometry of the graph G. 
Proposition 6.6 If D = 0, then system (|6| does not admit non-constant periodic solutions. 
Proof : 

Every point of a periodic solution is a limit point and, according to our preceding results (corollary 6.4), 
it is a critical point. Unicity of a solution provides a contradiction. □ 

Proposition 6.7 If D ^0, then the set of critical points of system ([6| is non- countable. 

Proof : 

Let X]fe=i'^fc(0) — G > 0. We know that the corresponding solution has to remain in the hyperplane 
(H) : X]fc=i ■'^k — G. Since the path is bounded it admits at least one limit point and, according to our 



preceding results (corollary 6.4), the latter is a critical point belonging to (H). Consequently, for every 



positive value of C, we obtain distinct critical points. □ 
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7 Stability of pure transport processes 



7.1 The irreducible case 

We consider pure transport processes (i.e. D = 0) on general graphs. We first discuss the stability of the 
special class of critical points a > solving equations of the form Fa = cl. Without loss of generality, 
we set T ~ 1. For such a, iVj(a) = N — k + c, and therefore, when the graph is regular, one obtains for 
example the uniform solution a = (p) = {p, . . . , p). When G is the complete graph of L nodes, where 
every pair of nodes i j are nearest neighbours, a simple computation shows that the Jacobian d/((p)) 
associated with ^ and evaluated at the uniform configuration (p), is given by 

dMip)) ^ dMip)) ^ ^ p' 

Consequently, d/((p)) is a symmetric generator, and thus admits only non-positive real eigenvalues. Fhe 
uniform configuration is then stable for the complete graph. 

Proposition 7.1 Let a > be such that Fa = cl, for some positive constant c > 0. According to Lemma 
a is a critical point, with Ni(a) = N = c + k. Assume that D — and set T — 1. The Jacobian 



5.3 



d/(a) — (dfi/daj) evaluated at a is then given by 

df{a)^^d{a)r{c id-d(a)F 
where d{a) is the diagonal matrix of diagonal given by a, and where F is the adjacency of the graph. 
The proof of Proposition |7.1| is given in Section |8] 

We now characterize the set of stable configurations using the spectral gap of the matrix P{a) = Td{a)/c. 
Let P be a stochastic matrix associated with a Markov chain on the state space A. We assume that P 
is reversible with invariant probability measure tt. Let A be the matrix defined by Aij = -KiPij = T^jPji, 
i ^ j. The eigenvalues of P are real, given by— 1 < < ■ ■ ■ ^2 < Pi — 1, and the spectral gap is 
given by C = 1 — /?2- Let L = id — P be the associated Laplace operator, of eigenvalues = 1 — (3^, 
fc = 1, • • • , L. Then (see e.g. ( Diaconis and Stroock,||1991 )) 



C = X2 = infj ^^^'^' : (j) is nonconstant}, (20) 
Var^(0) 



where 



is the Dirichlet form associated with L, and where Var^(0) is the variance of the random variable (j) with 
respect to the invariant probability measure tt. One can check that 



2 



We can also reformulate the above variational problem in a different way: set (0)^ = X)ieA Then 
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Lemma 7.2 Let G be a connected graph of adjacency matrix T, and let a > satisfy Ta = c 1 for some 
c > 0. The matrix P{a) defined by 

P{a) = ^rd(a), (22) 

c 

is stochastic, irreducible, reversible, of invariant measure 7r'(a) given by T:'{a)i = ai/{pL), and with a 
real spectrum —1 < /3a < /3a-i ^ • • • < /32 < /3i = 1. Let C{a) be the spectral gap of P{a), defined by 
C{a) = 1 — j32. a is stable if and only if C{a) > 1. Moreover, the spectral gap is given by 

where 6 = pL/c > 1, and where the infimum is taken over all nonconstant functions (p. 
Proof : 

The matrix is stochastic since by assumption Fa = cl. 

Let 7r'(a) — {jft) , ■ Then P{a) is reversible of invariant measure given by tt'. Notice next that 

Aij = Ti'{a)iP(a)ij = (57ij7r'(a)i7r'(a)j, 

where we recall that 7^ S {0, 1} is the entry of the adjacency matrix F. Hence, using the variational 
characterization of the spectral gap given in ( 20 1 , 

- Var^,(„)(</.) E.,,(^(i)-</'«)'T'(a),^'(a), " ' 

when (j) is non-constant. The configuration a is stable if and only if the eigenvalues of the Jacobian 
matrix df{a) given in the proposition 7.1 are all non-positive. The adjacency matrix F is symmetric, so 
that (Fc?(a))* = d{a)T. It follows that the eigenvalues f3i of d{a)T are equal to cPi, i = 1, ■ ■ ■ , L. The 
eigenvalues of N'^df{a) are given by jii{c — Pi) = /3i(l — /3i)c. Hence, a is stable if and only if (^2 < 0, 
that is if and only if C > 1. □ 

Corollary 7.3 Let G be a connected graph of adjacency matrix T, and let a > satisfy Ta — c 1 for 
some c > 0. For i Cz A, let Vi = {j G A; j ^ i} be the neighbourhood of i. Assume that there exist 
elements zq, ii, 12 o,nd i^ of A such that 

e i2 e V,, \ \ {io}, is e V,, \ \ {zo}. (23) 

Then a is unstable. 

Example 7.4 When G is a sub-graph of a two-dimensional grid, a solution to the linear system Ta = cl 
can possibly to be stable only when G belongs to the list given in Figure [1| which consists in the square, 
the star, and all the various parts of the star. 

Proof : 



We use Lemma 7.2 to express the spectral gap of P{a) as 

C{a) = 5 mf 

^ .^^ E» <l>{i?^''{c^)^JL - lij(f>ii)<j>{j)n'{a)y{a)j 



S inf 



0(z)^7r'(a),<5-i - E.,, l^MmJy{ah^'{a), 
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We will prove that C{a) < 1 by choosing a test function (f) satisfying {(t>).K'{a) = for which 

Ei<^«^T'(a), < ' 

which is equivalent to require that 

^lij4>{i)4>{iW{a),Ti'{a)j > 0. 

We set — 0, Vj G ViQ. For j e A \ Vig \ {io}, we choose to be arbitrary but positive. For j = io, 
we choose (/)(io) so that 

Consequently {<i))-n'(a) = and li3^i'^)'PU)aia] > 0. □ 

Corollary |7.3| provides a simple condition ensuring the non-stability of configurations a satisfying Fa — cl. 
We next consider the reducible case where = for i g / C A. Set J — A \ I, and let {71, • • • , jp} be 
the collection of sub-graphs of G induced by the nodes of J, of node set J^^ and of adjacency matrices 
^jpi P = li ■ ■ ■ J P- We again assume that F^^aj^^ = c^pl for some c^^ > 0. 



7.2 The reducible case 

We consider the stability of critical points a such that ai = 0, for i e / C A with / 7^ 0. 

Proposition 7.5 Assume that D ~ and set T = 1. Let a be a critical point of Q such that ai ~ 
for i Cz I . Let {7i,--- ,7p} he the collection of sub- graphs of G obtained by deleting the nodes of I , of 
adjacency matrices F^^, p — 1, - ■ ■ ,P. The critical points a are obtained by solving linear systems of the 
form F^pttl^p = c^pll^j, for s ome > (see Section |5. The spectrum of the Jacobian evaluated at a 
is given by 

spec(d/(a)) = U spec (d/|,^(a|,J) U jj] ^ - e /| (24) 

The proof of Proposition |7.5| is given in Section |8] 

Proposition 7.5 shows that such configurations are stable when 1) each a\^^ is stable and 2) when 
Efc-j afc/7Vfc(a) - (iVj(a) - K)/Ni{a) < 0, i £ I. Here, if fc - i, i € /, fc G A \ /, Nk{a) is given by 
the constant k + Cj^ when k € Jj^ . To go further, we need the following 

Definition 7.6 Let J C A. The outer boundary of J, denoted by dJ, is the subset of A given by 

dJ = {jeA\J; J}. 



7.3 Example: the rectangular grid 

We now illustrate the various stable patches we can form by using the building blocks, as given in Figures 
[2] and [3] It is easy to provide examples of unstable configurations when the outer boundary of some 
component 7 is such that 

d(^dJj^ n Jy 7^ 0, for some component 7' 7^ 7, (25) 

as illustrated in Figure |5ja) . 
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i ^ I 



■ • ■ ■ • ■ 

■ • ■ ■ o ■ 



• • • ■ • 

(a) Unstable configuration 



• • • • ■ • 

(b) Stable configuration 



Figure 5: (a) One can check in this example that (251 implies the non-stability of the configuration for 
well chosen parameters. Red dots indicates cells with 7^ 0. (b) One can check in this example that 



(26 1 is satisfied, ensuring the stability of the configuration. 



Next, the reader can verify, using Proposition |7.5| that any patch composed of building blocks disposed 
in such a way that 

d(dJp u Jp) n ( Up,^p Jp.^ = 0, Vp = 1, • • • , p, (26) 

is stable. Figure Islb) exhibits a typical example of a stable configuration in this setting. 



7.4 Example: the pure transport process on the circle 



We here assume that D = Q and T = 1. Corollary 7.3 yields the instability of uniform solution {p) — 
(p, . . . , p) when the length L of the cycle is larger than 4. The adjacency matrix of the circle is circulant, 
with eigenvalues given by 



= 2 cos 27r 



L 



The determinant of F vanishes if and only if there exists j € {1, L} such that /ij =0, that is if 



cos 27r 



^ 2ttj = ^ + /ctt for fc e N, 

or equivalently if there is a fc € N such that j — ^ + £ N. Hence, the determinant of F vanishes if 
and only if L is a multiple of 4. In this case, the set M^, of critical values a (that is satisfying Fa = cl) 
such that > 0, Vi G A, is such that 

= c — fli , 04 = c — a2 , 05 — 01,06 — 02,(17 — c — ai, ... 
with oi 7^ 7^ 02. Recalling that we impose the following normalization X^i^i '^i ~ P^i obtain 

L 



"^Ui = pL <^ 2c— = pL ^ c = 2p 



i=l 

The set of critical values Mc is then composed of configurations of the form 

a — (oi, 02, 2p — ai,2p — 02, oi, 02, 2p — ai,2p — 02, oi, 02, 2p — ai,2p — 02) 



with (01,02) G (0,2p) X (0,2p). Corollary 7.3 then implies that this set contains only unstable points 
when L > 4. For L — 4, the critical point a = (oi, 02, 2p — oi, 2p — 02) is stable since the eigenvalues of 
the Jacobian matrix are such that 

2c 



Ai = A2 = A3 = and A4 



(k + c)2 
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We can summarize these results in the following corollary: 



Corollary 7.7 Assume that the nodes are arranged on a circle of size L. The set Mc of critical values 
a > such that f{a) = contains only the uniform configuration (p, ...,/o) if L is not a multiple of 4- In 
the case where L — An, for some n G N with n > 1, is given by 

Mc = {(ai, 02, -fli + 2p, -02 + 2/9, 01,02, -fli + 2p, -02 + 2p, • • • ); ak € (0, 2p), k = I, 2}. 

Any element of Mc is unstable except for L = A. 



The set M*°* of all critical points is obtained by decomposing the circle into sub-graph 7 such a|-y > 
and by solving the system 



I^.-y Ql I ^ Cf-y 1 I <-y ^ 



for these sub-graphs. We can prove that this system has positive solution if and only |/| < 4 (|/| := 
length of the path), because for |/| > 4, we see that 04 = (which is in contradiction with the hypothesis). 
When |/| = 3, the critical points take the form a\.y — [zi,c^^c^ — zi), with zi G (0, c^) and when |/| — 2, 
oL = (c-y,c-y). In these two cases, the critical points are stable as the Lyapunov function H defined in 



(18 1 takes its minimal value H{a) = —{k+ ^)pL. The global minimum of H is obtained by adapting 
the result of Motzkin and Straus, (1965), see Remark 6.5 Finally, if |/| = 1, we have a\j — (c-y); H is 



maximal and hence a is unstable. 

The set M*°* of critical points is then obtained by taking the direct product of the sets of critical values 
associated with the paths 7. For example, if L is a multiple of 4, the subset of M*°* defined by 

Mc = {(oi,02, -oi -I- 2p, -02 4- 2p, oi, 02, -oi + 2p, -02 -f- 2p, • • • ); Oi = 0, 02 € (0, 2p)}, 

is composed of critical values which are stable since 

L —2(? L 

A = with multiplicity 3— and A = ^ with multiplicity — 

4 [k + c)^ 4 



7.5 An explicit computation when D = Q on the circle 

As we have seen, when |/| =3, the stable configurations are given by triplets of the form (zi, c-y, - 
where zi is such that zi G (0,c^), for some positive constant > 0. 

Consider a path composed of five cells i ~ 1, i, i + 1, i + 2 and i + 3 such that a^-i = 0^+3 = 0, so 
the dynamical system (|2| associated with these cells becomes 



doj 

"d^ 
dai+2 

dt 
doj+i 

dt 



Oi+iO; 



K + Oj + O, 



1+2 



OiOi+l 
1^ + O-i+l 



K + ai + 0^+2 K + Oj+l 

K + Oi+1 K + Oj+l K + Oj + Oi+2 K + Oi + 0,^+2 



that 

(27) 
(28) 
(29) 



Dividing (27 1 by (28 1 yields that 



dt 
dt 



ai+2 



Thus there is a positive constant c > such that 



Oi+2 = COi 



(30) 
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Plugging this identity in ( 29 1 , one obtains 



dt 



(1 + c)a^ai+i{ 



1 



and finally 



dt 

dt 



-(1 + c). 



Hence there exists a constant d such that ai+i = d — (1 + c)ai. NorniaHzing the total mass in such a way 
that Ui + ai+i + ai^2 = 3p, one gets that 3p — d and 



a,+i = 3p - (1 + c)ai. 



(31) 



Plugging (30) and (31 1 in equation (27 1 yields the differential equation 



doi 
"dT 



aj(3p - (1 + c)ai){3p - 2(1 + c)a^) 



(k + ai{l + c))(3p + K — (1 + c)aj) 
Setting w = (1 + c)aj, one gets the o.d.e. 



du 
dt 



u{3p — u){3p — 2u) 



(k + u)(3p + K - 
Solving by partial fractions expansions, one obtains 

Skp + k^ / X wo 9/92 + 4(k2- 



-(ln(M) + ln(3p- u)) 



3/9k) 



ln(3p — 2u) — t + a, 



9p2 ^ ^ ^ ^ 18p2 
for some constant a. Clearly one must have u < 3p/2. 
Lemma 7.8 As t ^ oo, u{t) = (1 + c)a^{t) — > ^. 
Proof : 

The preceding considerations show that we have to consider only initial conditions of the form < m(0) < 
3p. Clearly 0, ^ and 3p are critical points of our equation. 

We can easily find a compact interval / whose interior contains J = [0, 2>p] and so that f'{u) is continuous 
and thus bounded over /. As a consequence / satisfies a Lipschitz-condition over I. According to 
the general theory, for any initial condition u(0) g J our equation admits a unique solution defined 
over a maximal interval Im- If w(0) — 0, then m = is the corresponding solution. If u(0) e]0, then 
w(0) > 0. Due to unicity, the solution can not reach a critical point in a finite time and thus the boundary 
of ]0, Moreover the solution is obviously bounded entailing /,„ — [0, +co[. For the preceding reasons 
the derivative of u{t) is never and thus always positive since ?i(0) > 0. Thus u{t) increases to ^ as 



t — >■ +O0. The same reasoning shows that u{t) decreases to ^ as t ^ +oo for m(0) G]-^,3p[. Finally if 
u(0) = 3/9, then u = 3p. □ 



Furthermore, (30 1 yields 



H+2 



(0) 



a,(0) 



As (1 + c)ai 



tends to as time goes to infinity. Lemma 7.8 yields that = 3p/2, and 



ai{t) 



3p 
2(1 + c) 



1 + c 
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as i CXI. (30 1 and (31 1 show that 



a,;+i = 3/3 — (1 + c)a,: — > — Cy and 0^+2 — cui — > - — ^ — Cj 



2 ' ' 1 + c ' ' 1 + c 

In summary, one obtains that an orbit defined by initial conditions of the form 

(ai_i(0),ai(0),ai+i(0), 0^+2(0), a.j+3(0)) with ai_i(0) = aj+3(0) = 

converges to the critical point (zi,c-y,c-y — zi), with zi = = ^ and c = °a^(o)'^ ■ Finally, if the 

system starts from a symmetric initial state 0^(0) = 0^+2(0), the constant c is egal to 1 and the system 
tends to (0, ^, ^,0) as t ^ 00. 

8 Appendix 

8.1 Proof of Theorem l3Jl 

First, we easily check that the system a = f{a) is conservative, i.e. 

L L 

Vi eM>o, ^a.(t) = ^a,(0). 

i i 

In the following, we use the notation a instead of The latter is equivalent to 

L L 

i i 

In fact, one can write 

i i k^i i k'^i 

= D Y,{d^a^ - d,a,) + 2T^ ^^((7V, - Nk) - [Nk - N,)) = 0, 

where Nk = k + ^k, and where di is the degree of i (that is the number of neighbours of i). 

Next, system Q can be written as 



a^. (32) 



Let a, a solution of (^32| with ci(0) G II^>o- 

We say that the function / : M_|_ — > M is instantaneously positive (i.p.) if there exists (5 > so that / 
is strictly positive over (0,(5). If /(O) > and / is continuous to the right at 0, then / is i.p.. It is also 
clear that if / admits a strictly positive right-hand derivative at 0, then it is i.p.. 

Let U be the open set U — {x — {xi,X2, ■■■,xl) e M^;— ^ < Xi}. Since the right-hand member of 



(32 1 is continous over U, the general theory of o.d.e.'s provides the existence of a solution defined over 
a maximal interval G J+ C M+ for any initial condition a(0) G U. Moreover, the solution is unique 
because the right-hand member of (p2|) locally lipschitzian. Set for convenience 



hi{t) = Z)^^a^,(t) and 

k^i 

ak{t) ak{t) D 

T 



(f\ = r V ( °fcft) Qfc(^) 
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The variation of constants formula allows us to write , Vt G J"*", 



Since a,;(0) > 0, the first term in ( |33[ ) is non-negative. Moreover if ak{t) is i.p. for some k ^ i, then 
according to (33 1, the same property holds for ai{t). In particular, if afe(O) > for some k ^ i, then by 
continuity ak(t) is i.p. and thus also ai{t). 

The case D > 0: 

Clearly, if a(0) = 0, then the unique solution is identically 0. Otherwise, there exists 1 < io < L with 
Ojp(O) > and Vj io,aj(t) is i.p.. Since our graph is supposed to be connected, every i admits a 
neighbor k ^ i with afe(t) i.p.. Hence, ai{t) is i.p. Vi,! < i < L. 

The preceding arguments show that for any initial condition a(0) G M>q C U, all components of the 
solution of (32 1 are i.p.. Let us suppose that one of them admits the value in J"'"\{0}. Since all 
components are continuous and their number is finite, there exists a first time tg > for which at least 
one component aig{to) = and all of them are strictly positive over (0,io)- According to (33 1, we have 

Jo 

Clearly hi{t) > over J+\{0} and since the first term is non-negative, we conclude to ai„(to) > 0, a 
contradiction. Therefore all ai{t) are strictly positive over J"'"\{0}. 

The case D = 0: 

If ai(0) — 0, the homogeneous equation for a.i(t) admits only the zero solution, and we remove the 
related zth component from (32 1. Otherwise aj(0) > and, by continuity, ai{t) is i.p.. In that case 
a,it) = a,(0)e-/'o > o over J+. 



In both cases the solution of (32) have strictly positive components over 
\/t G J"*" we have: 

^ a,{t)= a^(0)- 

l<i<L l<i<L 



We also proved that 



As a consequence the solution of ( 32 ) is bounded and thus the unique solution of our problem is defined 
over J+ = [0, +oo). 



8.2 Proof of Proposition 7.1 



We first give the Jacobian, for general a. We have 

dfi{a) _ at _ fli 



da, N, 



(where the last term is due to the triangles in the graph) when j ^ i, that is, i and j are nearest 
neighbours. When i — j, one gets 

dfi{a) _ _ sr^ ak_ I]fc~i "-k 



k^i 



ak_ 
Nk 



N, 



(35) 



The remaining non-vanishing partial derivatives correspond to nodes j located at distance 2 of i in the 
graph, that is, to nodes j such that j ^ k for some k ^ i, j i but i j- Then 



9a, 



E 



j^k, k^i 



aittk 



(36) 



22 



When Ni = N, Vi, these expressions simplify to 



dfi{a) sr-^ UiGk N~K a 

If 



daj ~^ ~ 7V2 ^2 E "fc- 



and 

dfi{a) _ \ - giflfc _ \ - 

if j k for some k ^ i, j i but i 7^ j. 

Consider the sub-matrix i given by L = {dfi{a)/daj)j^i. Let rf(a) be the diagonal matrix of diagonal 
given by a. The perturbation associated with the triangles contained in the graph is represented by the 
term —j^ Y^kf-^i kr^j '^k in for j i, and the related matrix is given by 



J2 «fe = ^likaklkj J lij 



id{a)Td{a)T - diag(d(a)rd(a)r)) . .J 7,,- 
1 TV - K \ 



The matrix L is now given by 

d{a) 



L = 



^2 - - ^'i) - ^(d(a)rd(a)r - (iV - «)(i(a)) o r, 



where o represents the Hadamard product, i.e. the multiplication component by component. 
Likewise, the perturbation of L by ( ) can be written as 

V 3 / ir^k,k^j,i'/'j,i¥^j 



Ar2 



k^i^k^j, . 



lij = ^ X] 'TifcCtfeTfej^ (1 - lij - idij) 



1 N — K 

= ( -Jpid{a)'^dia)r),j + d{a)ij ) (1 - 7,^ - 



The related Jacobian is thus given by L + ( ^g^"'^ ) , that is 

V ^ / ir^k,k^j,i'7^j,i¥^j 

^/W = ^iN-K){T-id)-^id{a)Td{a)T-{N-K)d{a))oT 



id{a)rdia)r -{N- K)d{a)) o{l-T-id) 
1 



1 

^{N- k){T - id) - ^{d{a)Td{a)T - {N - K)d{a)) o (1 - id) 
^{N- «)(r - id) - ^{d{a)Td{a)T - (TV - K)d{a)), 
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where 1 is the matrix composed only of ones. The last equality is a consequence of the fact that the 
diagonal of d{a)Td{a)r — (N — K)d{a) vanishes. Hence, 

^/(«) = ^((^ - - d{a)r) ^ - d(a)r), 

proving the result. 



8.3 Proof of Proposition 7.5 



Set / = {i e A : a,; = 0}, and consider the sub-graphs 7^ of G induced by the nodes of J = A \ /, with 
7p — {Ap,Ep), 1 < p < P. The related critical points a are such that the restrictions a\^^ satisfy the 
linear systems T^^al^^ ^ Cj^llj^. Set Nj^ = Cj^ + k. 

(34 1 - (36 1 permit to compute the entries of the Jacobian matrix, by first looking at the diagonal entries: 
When i G Ap, one has 

da, - Nl ' 

providing the diagonal entry of the Jacobian of f\'yj,{a\j^). When i ^ Ap, a similar computation yields 

dfija) _sr^ 0^ _ - K 
da, ~f^.Nk N, - 

We then compute the entries for j ^ i: 

dfija) _ - K ai_ _ - ^ _ aka-i 

for i,j e Ap and 1 < p < P, which corresponds to the (i, j) entry of the Jacobian of f\~f^{a\^^). Likewise, 

J J rp 7p k~'i,k~j, 'p ■' <p <p fe~j,fc~j, 



when i G Ap for some p and j ^ Ap. Finally, 

a/, (a) 



9aj 



0, 



when i, j ^ UpAp, or equivalently when both i and j belongs to /. 

We next consider (i, j) entries where j is at a distance 2 of i in the graph G, that is when j is such that 
j ^ k for some k ^ i, j ^ i and j 'j^ i. One obtains that 



when i, j, fc € Ap, which is the {i.f) entry of the Jacobian of ]\^^(a\^^ . 
Likewise, 

da ~ ^ iV2 ' 

j~k,k~i '>P 

when i, A: € Ap, j ^ Ap j e I). 
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Next, 

ma) ^0 

daj 

when i or fc ^ Ap, Vj G A. 

Permuting conveniently the indices, the Jacobian d/(a) can be written as 

dp) (37) 

where d„ is a diagonal matrix n x n with entries given by Ai := X]fc~j nZ ~ ^n'-'" ' "^'^'^ i & I, and hence 
d/'^ is a block diagonal matrix, each block being equal to the Jacobian of / restricted on each sub-graph 
7p. The permutation allows us to group all indices i e / in the same block, and all indices related to 
the sub-graphs 7p are also arranged together. It follows that the eigenvalues of d/(a) are given by the 
diagonal entries (Ai)ig/, and by the eigenvalues of all Jacobian matrices. 
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